#!/usr/bin/env python
from math import sqrt, hypot, acos

class Point:
    def __init__(self, x, y):
        self.x = x
        self.y = y

class Line:
    def __init__(self, x0, y0, x1, y1):
        self.x0 = x0
        self.y0 = y0
        self.x1 = x1
        self.y1 = y1
        self.dx = x1 - x0
        self.dy = y1 - y0
    
    def collide_line(self, x0, y0, x1, y1):
        t = (x1 - x0) * (self.y0 - y0) - (y1 - y0) * (self.x0 - x0) / \
            (y1 - y0) * (self.x1 - self.x0) - (x1 - x0) * (self.y1 - self.y0)
        return self.x0 + self.dx * t, self.y0 + self.dy * t, t

class Curve:
    def __init__(self, x0, y0, x1, y1, x2, y2):
        self.x0 = x0
        self.y0 = y0
        self.x1 = x1
        self.y1 = y1
        self.x2 = x2
        self.y2 = y2
        
        self.length = self.getLength()
        
        left = min(x0, x1, x2)
        top = min(y0, y1, y2)
        right = max(x0, x1, x2)
        bottom = max(y0, y1, y2)
    
    def getLength(self, segments=32.0):
        """The length of a curve and the other variables calculated herein
        are extremely useful to calculate the projection of points on this
        curve as precisely and efficiently as possible.
        """
        self.length = 0.0
        
        maxLength = 0.0
        minLength = 255.0 #Won't ever be this high, really.
        
        p1 = self.p1.x, self.p1.y
        for i in xrange(1, int(segments)):
            p2 = bezier(self.p1.x, self.p1.y,
                        self.p2.x, self.p2.y,
                        self.p3.x, self.p3.y, segments / float(i))
            
            segmentLength = hypot(p1[0] - p2[0], p1[1] - p2[1])
            
            maxLength = max(maxLength, segmentLength)
            minLength = min(minLength, segmentLength)
            
            self.length += segmentLength
            
            p1 = p2 #Prepare for the next iteration
        
        self.maxVariation = maxLength / minLength
        
        return self.length
    
    def makeImage(self, segments=32):
        """For optimization purposes, curves are drawn once on a surface which
        is then blitted as needed.
        """
        image = pygame.Surface(self.AABB.size, pygame.SRCALPHA, 32)
        image.fill((0,0,0,0))
        image.convert_alpha()
        
        p1 = self.p1.x - self.AABB.x, self.p1.y - self.AABB.y
        for i in xrange(1, int(segments)):
            p2 = bezier(self.p1.x - self.AABB.x, self.p1.y - self.AABB.y,
                        self.p2.x - self.AABB.x, self.p2.y - self.AABB.y,
                        self.p3.x - self.AABB.x, self.p3.y - self.AABB.y,
                        segments / float(i))
            
            pygame.draw.aaline(image, (0,0,0), p1, p2)
            
            p1 = p2 #Prepare for the next iteration
        
        return image
    
    def getPoint(self, t):
        return bezier(self.p1.x, self.p1.y,
                      self.p2.x, self.p2.y,
                      self.p3.x, self.p3.y, t)

#Vector operations
def sum(v1, v2):
    return v1[0] + v2[0], v1[1] + v2[1]
def dot(v1, v2):
    return v1[0]*v2[0] + v1[1]*v2[1]

#Line operations
def lineEquation(Ax, Ay, Bx, By):
    """Implicit line equation, returns A B and C
    """
    return By - Ay, Bx - Ax, Ay * Bx - Ax * Bx

def lineIntersection(((x1, y1), (x2, y2)), ((x3, y3), (x4, y4))):
    """Intersection point of two lines
    """
    delta = float((y3-y4)*(x2-x1) - (x3-x4)*(y2-y1))
    if delta == 0:
        return None
    t1 = ((x3-x1) * (y3-y4) - (y3-y1) * (x3-x4)) / delta
    t2 = ((x2-x1) * (y3-y1) - (y2-y1) * (x3-x1)) / delta
    x = (x2-x1) * t1 + x1
    y = (y2-y1) * t1 + y1
    segments = False
    if (t1 >= 0.0) and (t1 <= 1.0) and (t2 >= 0.0) and (t2 <= 1.0):
        segments = True
    return x, y, segments

def linesCollide(((x1, y1), (x2, y2)), ((x3, y3), (x4, y4))):
    """Just returns whether two lines collide or not
    """
    delta = (y3 - y4) * (x2 - x1) - (x3 - x4) * (y2 - y1)
    if delta != 0:
        t1 = ((x3-x1) * (y3-y4) - (y3-y1) * (x3-x4)) / delta
        t2 = ((x2-x1) * (y3-y1) - (y2-y1) * (x3-x1)) / delta
        if 0.0 <= t1 <= 1.0 and 0.0 <= t2 <= 1.0:
            return True
    return False

def pointOnLine(Px, Py, Ax, Ay, Bx, By):
    """Exactly like the above, but it also returns d
    """
    Vx = Px - Ax
    Vy = Py - Ay
    Cx = Bx - Ax
    Cy = By - Ay
    d = (Vx * Cx + Vy * Cy) / (Cx * Cx + Cy * Cy)
    return Ax + Cx * d, Ay + Cy * d, d

def pointOnSegment(Px, Py, Ax, Ay, Bx, By):
    """Exactly like the above, but it also returns d
    """
    Vx = Px - Ax
    Vy = Py - Ay
    Cx = Bx - Ax
    Cy = By - Ay
    d = (Vx * Cx + Vy * Cy) / (Cx * Cx + Cy * Cy)
    if d > 1.0: d = 1.0
    if d < 0.0: d = 0.0
    return Ax + Cx * d, Ay + Cy * d, d

#Curve operations
def bezier(Ax, Ay, Bx, By, Cx, Cy, t):
    """Calculate a given point on a quadratic bezier curve.
    
    @param t: a float between 0 and 1
    """
    x = (1 - t)**2 * Ax + 2 * t * (1 - t) * Bx + t**2 * Cx
    y = (1 - t)**2 * Ay + 2 * t * (1 - t) * By + t**2 * Cy
    return x, y

def lineCurveIntersection(Ax, Ay, Bx, By, Cx, Cy, Dx, Dy, Ex, Ey):
    """A and B are the points of the line, C D and E are those of the curve"""
    a = (By - Ax) * (Cx - 2*Dx + Ex) + (Ax - Bx) * (Cy - 2*Dy + Ey)
    b = (By - Ax) * (-2*Cx + 2*Dx) + (Ax - Bx) * (-2*Cy + 2*Dy)
    c = (By - Ax) * (Cx - Ax) + (Ax - Bx) * (Cy - Ay)
    squareDelta = b**2 - 4 * a * c
    if squareDelta > 0:
        delta = sqrt(squareDelta)
        t1 = (-b + delta) / (2.0 * a)
        t2 =  (-b - delta) / (2.0 * a)
        p1 = bezier(Cx, Cy, Dx, Dy, Ex, Ey, t1)
        p2 = bezier(Cx, Cy, Dx, Dy, Ex, Ey, t2)
        return p1, p2, t1, t2
    elif squareDelta == 0:
        delta = sqrt(squareDelta)
        t = (-b + delta) / (2.0 * a)
        p = bezier(Cx, Cy, Dx, Dy, Ex, Ey, t)
        return p, p, t, t
    else:
        return None

def pointOnCurve(Px, Py, Ax, Ay, Bx, By, Cx, Cy, approximations=4):
    lowerLimit = 0.0
    upperLimit = 1.0
    x0, y0 = bezier(Ax, Ay, Bx, By, Cx, Cy, lowerLimit)
    x1, y1 = bezier(Ax, Ay, Bx, By, Cx, Cy, upperLimit)
    for i in range(approximations):
        if hypot(Px - x0, Py - y0) < hypot(Px - x1, Py - y1):
            upperLimit = (upperLimit + lowerLimit) / 2.0
            x1, y1 = bezier(Ax, Ay, Bx, By, Cx, Cy, upperLimit)
        else:
            lowerLimit = (upperLimit + lowerLimit) / 2.0
            x0, y0 = bezier(Ax, Ay, Bx, By, Cx, Cy, lowerLimit)
    return (x0 + x1) / 2.0, (y0 + y1) / 2.0